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Abstract 



Accelerator design usually requires that the betatron tune is kept away from key values which would 
allow for the growth of transverse particle oscillations due to resonant conditions driven by lattice 
imperfections (the consequence of which being that particles with increasing transverse amplitude 
will eventually collide with the accelerator wall). The non-scaling fixed field alternating gradient (NS- 
FFAG) accelerator allows the betatron tune to vary with particle momentum and in doing so to cross 
conventionally forbidden values. 

This dissertation project looks at the approximation of the proof of principle linear NS-FFAG 
accelerator EMMA (Electron Machine with Many Applications) to a constant gradient accelerator, 
with the aim of simulating the transverse amplitude growth that occurs as a result of crossing 
integer betatron tune values and transverse misalignment in quadrupole positions. 

The influence of the magnitude of the quadrupole alignment errors and the rate of acceleration are 
studied, and the results from this approximate model are compared to those obtained using a 'full 
physics' model 1 . 

Results 

• The constant gradient approximated model was found to be able to demonstrate the same 
features of integer tune crossing as a 'full physics' model. 

• The amplitude growth as a function of acceleration rate was found to have a gradient of 
(38 + l)m _1 /2 (3% error). 

• The success of a linear NS-FFAG in being able to accelerate particles to extraction energy, 
was found to be dependent upon the rate of acceleration. 



Machida, S. (2008) Resonance Crossing and Dynamic Aperture in Non-scaling Fixed Field Alternating Gradient 
Accelerators 

'full physics 7 is here defined as model which does not make the approximation of alternating gradient focusing to being 
constant gradient. 
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Section I - Introduction 



The Non-Scaling Fixed Field Alternating Gradient (NS-FFAG) accelerator is a novel accelerator design 
upon which the British Accelerator and Radiation Oncology Consortium (BASROC) intend to base a 
proton therapy facility. The NS-FFAG has the potential to give both economical and practical 
advantages over cyclotron and synchrotron alternatives when applied to the field of oncology. 



Radiation Therapy for Treating Cancer 



The first page of the British Accelerator Science and Radiation Oncology Consortiums (BASROC) 
website claims that one third of all UK residents will, at some time, suffer from cancer, and that it 
can be expected that the treatment for half of these people will involve radiation therapy 2 . 

Radiation therapy (or radiotherapy) utilises ionizing radiation for the aim of destroying cancerous 
cells. Tumours made up of cancerous cells are found in the midst of healthy, functioning tissue, and 
so there is the constraint upon the application of radiotherapy that healthy cells must be preserved 
to the greatest possible degree. In fact, an effective treatment is often a compromise between the 
probability of eradicating the cancer and the probability of complications due to normal tissue 
damage 3 . 

Teletherapy (meaning distant therapy) involves a beam of particles emitted from a source being 
directed at a tumour. The dose (which may be defined as the energy deposited per unit mass JKg' 1 ) 
given to a specific volume within patient is determined by the way in which the beams particles 
deposit energy along the beam path 4 . The characteristics of this deposition of energy are specific to 
the type of particle used (see graph 1). 



2 

BASROC. (2006) British Accelerator Science and Radiation Oncology Consortium. Available at 
http://www.basroc.org.uk/index.htm 
Fenwick, J. (2010) Biological Modelling. PHYS384:Radiation Therapy Applications Notes. 
Wu Chao, A. & Tigner, M. (1999) Handbook of accelerator physics and engineering. 
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Graph 1: Depth-dose curves for photons, electrons and protons. 



Currently electron and photon based treatments dominate the field of radiotherapy; however the 
use of other particles may offer distinct advantages. Below is a brief summary of the features seen 
for the particles detailed in the depth dose curves of graph 1: 



Electrons for clinical use may be accelerated, by linear accelerator, from low energies to energies 
suited to therapeutic application (typically within the range 4MeV to 20MeV). An electron beam will 
give a relatively uniform dose with depth until rapidly falling off (at approximately 2cm for the 4MeV 
beam), giving a close to ideal dose distribution for superficial tumours, with little damage caused 
beyond the depth of the tumour. 

An electron loses energy at a greater rate towards the end of its range, however the low electron 
mass means that it is readily scattered through large angles and that the tissue depth is not 
synonymous with the electron path length. It is for this reason that the uniform dose with depth is 
realised. 



5 Paul, H. (2009) Depth Dose Curves. Available at http://en.wikipedia.Org/wiki/File:Depth Dose Curves.jpg 



Electrons: 
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EJlotons: 

Therapeutic photon production may be based on the same linear accelerators used in generating 
therapeutic electrons, the difference being that the accelerated electrons are incident upon a high 
atomic number metal target (usually tungsten). 

Photons are classified by the potential through which the electrons used to generate them are 
accelerated, e.g. 6MV photon spectrum created through Bremstrahlung emissions as 6MeV 
electrons are incident upon a tungsten target. The offset in the maximum of the depth-dose curve 
from a Ocm depth (see graph 1) is dependent upon the energy of the beam (with high energy beams 
exhibiting greater offsets). This is a useful feature as it means limited dosage at the surface (skin, 
hence reducing patient discomfort), and a potentially higher dosage at the tumour depth. Beyond 
the maximum, the dose falls off exponentially with depth. 

Taking a small tumour, starting at a depth of 1.5cm and extending through to 2.5cm, as an example, 
it can be seen from graph 1 that treatment with a 6MV photon beam would result in a significant 
dose being delivered before and after the tumour. 

Protons: 

The rate at which a proton deposits energy in tissue increases rapidly shortly before it comes to rest, 
the relatively large proton mass (when compared to that of an electron) leads to little scattering 
within tissue, and so tissue depth and proton path length are similar. A beam of monoenergetic 
protons incident upon a patient therefore results in the feature marked as the Bragg peak in graph 1. 

The Bragg peak allows for a deep seated tumour to be treated whilst causing very little damage to 
surrounding healthy tissue. The rapid fall off of dose beyond the proton range makes proton therapy 
particularly well suited to tumours that are in close proximity to sensitive organs (such as the spinal 
cord). 

The relationship between a protons energy (E) and its range (R) in tissue may be described 

approximately by R t i S sue — [t^j > f rom this, it can be estimated that energies of 220MeV would be 

required to treat tumours at depths of up to 30cm. To date, the acceleration of protons for medical 
application has been carried out using circular (rather than linear) accelerators 6 . 

6 This has been attributed to our current technical abilities, and the relatively high cost-performance 
ratio that would be involved in accelerating protons to energies in the region of 200-250MeV using 
linear accelerator. Mayles, P., Nahum, A., & Rosenwald J.C. (2007) Handbook of radiotherapy 
physics: theory and practice. 
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Figure 1: Above left, a simulated distribution of dose to a tumour (circled by red line) and 
surrounding tissue when a photon beam treatment is applied (beams entering from multiple 
angles in order to spread the dose to normal tissue and reduce the risk of normal tissue 
complications). Above right, the application of a proton based treatment to the same tumour. 
Blue regions are 'cold' (receive very little dose), whilst colours towards red end of spectrum 
indicate increased dose. It can be seen that for the proton based therapy, the higher dose regions 
conform more tightly to the tumour volume, leading to an overall lower dose for normal tissue. 



7 Lemoigne, Y. (2009) Radiotherapy and Brachytherapy [electronic book]. 
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Existing Accelerator Based Proton Therapy Units 



1. The Douglas cyclotron located at Clatterbridge Hospital (Merseyside, UK) is able to deliver 
62 MeV protons for therapeutic purposes. A proton of energy 62 MeV corresponds to a 
range of approximately 30mm in tissue, which is adequate to treat tumours located within 
the eye (this is the Douglas cyclotrons current application, with the Bragg peak meaning that 
an ocular tumour may be treated effectively, whilst the risk of damage to the optic nerve is 
minimal). 

The output of the Douglas cyclotron is fixed at 62MeV which presents a problem; the Bragg 
peak does not necessarily occur at the most desirable depth (i.e. in the vicinity of the 
tumour). The solution comes in the form of perspex (a tissue equivalent material) blocks, 
which are introduced in to the beams path in order to ensure that the Bragg peak is found at 
a depth consistent with the tumour depth. 

Further to the above problem, the depths over which a tumour extends can often be of a 
greater range than the pristine Bragg peak (the Bragg peak due to monoenergetic particles) 
is able to cover, for this reason there is a need to 'spread' the Bragg peak: 



Normalised 

absorbed dose 

(%) Resultant spread peak 




Graph 2: (Normalised dose Vs Depth) In order to effectively treat a tumour it may be necessary to 

spread the pristine Bragg peak. 

This spreading of the Bragg peak is achieved through the use of a modulator, a perspex disc 
of stepped thickness which rotates through the beam line (such methods are referred to as 
passive scanning): 
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Figure 2: A perspex modulator. 

2. The Heidelberg Ion Therapy Centre (HIT) is able to offer proton energies within the range of 
48 to 221MeV, corresponding to a depth of 20 to 300mm in water 9 (as well as carbon ions 10 
of up to 430MeV), for therapeutic applications. The range of energies available at the HIT 
are arrived at through the use of a synchrotron. 

Cyclotron design sees the orbit of a particle growing from a very small radius to a much 
larger radius (proportional to the particles momentum, and inversely so to the magnetic 
field strength), whereas within a synchrotron particles are accelerated at a fixed radius of 
orbit. This radius is maintained by a magnetic field which increases throughout acceleration, 
giving an economical advantage to the synchrotron as a result of the relative sizes of the 
magnetic components. 

The synchrotron provides a flexible solution that allows the HIT to deliver a range of proton 
energies, and for dynamic scanning to be applied 11 . However, this flexibility comes at a cost; 
the necessity to ramp up the magnetic field strength in maintaining the orbit radius limits 
the synchrotron to providing a pulsed output of relatively (when compared to a cyclotron) 
low intensity. 



Bonnett, D.E.et al. (1993) The 62 MeV proton beam for the treatment of ocular melanoma at Clatterbrigde. 
Available from: http://bir.biriournals.org/cgi/content/abstract/66/790/907 

Mosthaf, J.M. et al. (2008) Overview of the communication structure of the hit accelerator control system. 
Available from: http://accelconf.web.cern.ch/Accelconf/pc08/papers/wep002.pdf 

10 

The application of which shall not be covered here. 
11 dynamic scanning is where the beam energy is varied so that perspex blocks no longer need to be 
introduced into the beams path. Dynamic scanning is preferred over passive scanning as it provides better 
conformation of dose to a tumour. 
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Figure 3: The Heidelberg ion therapy centre, from ion sources to treatment gantry. 
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Another Way: Fixed Field Alternating Gradient Accelerators (FFAG's) 



FFAG accelerators were first proposed in the 1950s. Unlike for a synchrotron, the magnetic field of a 
FFAG accelerator is not time dependent (hence 'fixed field'), but instead varies spatially with a 
higher field experienced at larger radii. This means that as a particle is accelerated within a FFAG, the 
radius of orbit varies, but to a much smaller degree than would be expected within a cyclotron. The 
alternating gradient refers to periodic areas of focusing and defocusing magnetic fields which are 
located around the accelerator ring to provide a net focusing effect (as is the case for a synchrotron). 
In terms of the accelerator aspects considered important for proton therapy, FFAG's exhibit the 
flexibility of the synchrotron and the beam intensity of the cyclotron. 



Parameter 


Synchrotron 


Cyclotron 


FFAG 


Beam Intensity 


Low 


Plenty 


Plenty 


Operation 


Not Easy 


Easy 


Easy 


Maintenance 


Normal 


Hard 


Normal 


Possibility of Carbon Ions 


Yes 


Expensive 


Yes 


Possibility of Variable Energies 


Yes 


No 


Yes 


Multiple Beam Extraction 


Difficult 


No 


Yes 



Table 1: Comparison of synchroton, cyclotron and FFAG features. 



FFAG's: Scaling Vs Non-Scaling 

The magnetic field of a scaling (S) FFAG increases logarithmically with orbit radius (described by 

fr\^ 12 13 

B = B 0 ( — ) , where k is the field index ), which has the advantage that the shape of a particles 

\T 0 J 

orbit is maintained throughout acceleration. The downside to an S-FFAG is that the magnetic 
components are relatively large, complicated and expensive. 

A FFAG design which allows for a more straightforward and compact magnetic structure is that of 
the non-scaling (NS) FFAG. With a magnetic field which increases linearly with orbit radius, NS- 
FFAG's allow the shape of a particle's orbit to vary with momentum. This variation in orbit shape 
could prove problematic. Resonant conditions which are usually avoided within the S-FFAG (such as 
integer betatron tune values (discussed in section II)) are met as a matter of course within the NS- 
FFAG. If resonant conditions are present for a sufficient amount of time, errors in the accelerators 
magnetic fields may stimulate increases in the transverse amplitudes of particles to the point that 
they are incident upon the accelerator wall. 



Edgecock, R. (2009) FF : AGs for Radiotherapy. IOP Mayneord Phillips Summer Schools Notes. 
Advantages of scaling orbit described under non-scaling FFAG. 
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The non time dependent magnetic field means that acceleration can occur rapidly within a FFAG (the 
rate at which a field can be ramped up is a limiting factor within a synchrotron). It is this high rate of 
acceleration that allows for resonant conditions within a NS-FFAG to be passed quickly without the 
particle beam being lost. 




Figure 4: The NS-FFAG allows the orbit shape to vary with particle momentum, blue and red 
coloured areas (labelled D and F) represent alternating gradient scheme 



BASROC has the aim of constructing a NS-FFAG type accelerator for the purpose of treating cancer 
patients through the application of proton therapy. This aim is anticipated to be realised in three 
steps 14 : 

1. EMMA (Electron Machine with Many Applications) is a proof-of-principle machine which will 
accelerate electrons from 10.5-20.5MeV. EMMA is currently under construction at the 
Daresbury Laboratory (April 2010), and once completed will be used to gain better 
understanding of NS-FFAG design and operation. 

2. PAMELA (Particle Accelerator with MEdicaL Applications), a 70-100MeV proton NS-FFAG 
prototype intended to show the practicality of the application of NS-FFAG's in the field of 
hadron therapy. 

3. The construction of a complete facility which is able to offer proton therapy. 



BASROC. (2006) British Accelerator Science and Radiation Oncology Consortium. Available at 
http://www.basroc.org.uk/vision.htm 
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Section II - Accelerator Concepts 



Beam Optics 



For a circular accelerator there is a clear need to affect a beam of particles in such a way that it 
maintains a well defined path around the circumference of the accelerator. For a beam of charged 
particles, their motion may be guided through the use of electromagnetic fields. The force applied to 
a charged particle within an electromagnetic field can be described by the Lorentz force: 

F = q(E + v x B) 

For an electron of energy 20 MeV (with velocity approximately equal to the speed of light) orbiting in 
a circular accelerator of radius 5m, the Lorentz force required to maintain the orbit is defined as 
being equal and opposite to the centrifugal force: 



P c = 



mv 



Which for the electron and accelerator described above, may be approximated to: 
„ ym 0 c 2 E 

r r ~ = — (a 20MeV electron is highly relativistic with a velocity very close to that of the 

-^/-* ^/-* 

speed of light) 

F c = 6.4 x 10~ 13 N 




Figure 5: Definition of accelerator coordinates. 



Equating with the Lorentz force gives: 



F L (x,z,s) = -(6.4 x 10 _13 ,0,0)JV 



Ft 

— = (4 x 10 6 , 0,0) = E + v x B 



Where q = -1.602 x 10 



-19 
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Which may be solved (non-exclusively) by: 

E = (4 x 10 6 , 0, 0) Vm' 1 and B = (0, 0, 0)mT 

or 

E = (0, 0, 0) Vm' 1 and B = (0, -13, 0)mT 
Where v = (0, 0, 2.998 x 10 8 )ras _1 

From the formula for the Lorentz force, it may be conceived that the orbit described may be 
maintained by having the electron travel through an electric field of magnitude 4 x 10 6 7m" 1 . 
However, the dielectric breakdown voltage of 3 x 10 6 7ra _1 for air 15 means that such an application 
is unfeasible. 

A magnetic field of 13 mT can yield the same force as the electric field of above, but without the 
complications of dielectric breakdown, it is for this reason that accelerator beams are predominantly 
manipulated using magnets. 16 

The relationship between Lorentz force and centrifugal force can be used to define a quantity known 
as the magnetic rigidity: 



mv s 2 



= qv s B 2 



V 

B z (x,z,s)r(x,z,s) = - 

q 

Where p is the particle momentum 

The Taylor series expansion of -fi z (x) gives: 

p 

-B z (x) =-[B z0 +^—x + -—- r x z + ••• 
p p\ ax 2 dx z 

It is the first two components of this expansion that will be covered within this section, with greatest 
precedence put on the quadrupole field. 

-B z0 = - - Dipole equation [2.1] 

= kx - Quadrupole equation [2.2] 

p dx 



15 



16 



Tarr, M. (2007) Dielectric Properties. Available at http://www.ami.ac.uk/courses/topics/0184 dp/index. html 
Wille, K. (2000) The Physics of Particle Accelerators an Introduction, p. 44 
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Dipole Magnets 

As already detailed, for a particle travelling in the direction of the s axis, a magnetic field parallel to 
the z axis of the accelerator will lead to a force being experienced parallel to the x axis. It is such a 
field which is applied in order to bend a particle around the circumference of an accelerator. Dipole 
magnets may be located at intervals around an accelerator in order to provide the field necessary to 
maintain a particles orbit. 

Quadrupole Magnets 

At any point within an accelerator a particle may demonstrate a trajectory which is divergent from 
that of the ideal particle, without intervention such a particle would be lost. The quadrupole, an 
arrangement of four magnetic poles which straddle the accelerator beam pipe, serves to focus 
divergent particles: 
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Figure 6: Sketch to show the arrangement of a quadrupole, with magnetic field and the associated 
force experienced by an electron (travelling into the page) marked. The orientation of this 
quadrupole will lead to focussing in the horizontal plane, and defocusing in the vertical plane. 

The symmetries of the quadrupole mean that a particle travelling through the origin of the x-z axis 
shown on the diagram will experience no force, whilst the magnetic field (and so force) increases 



linearly with respect to displacement from the origin for an off centred particle (i.e. 
constant). The direction of this force is dependent upon the orientation of the quadrupole. 



dB z 
dx 



The quadrupole shown in figure 6 will focus a particle displaced along the horizontal axis and 
defocus a particle displaced along the vertical axis (when considering an electron travelling into the 



Wilson, E. (2001) An introduction to particle accelerators. 
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page). A rotation of the quadrupole through 90° leads to focusing in the vertical plane and 
defocusing in the horizontal plane. 

The effect of a horizontally focusing quadrupole of length I upon a particle with initial horizontal 
displacement x 1 and divergence may be described by a transfer matrix: 



M QF = 



cos 



yfKl 



sin 



4ki 



—yfRsiny/Kl cos^l 



Where k = - 



q dB z 



p dx 



The displacement and divergence after the quadrupole (x 2 and x' 2 ) is given by: 



x^ 



The FODO Cell 

The FODO cell is a commonly utilised arrangement of quadrupoles which results in a net beam 
focusing effect in both the horizontal and vertical plane. The FODO cell consists of a focusing 
quadrupole, a non-focusing region (i.e. a drift space or dipole), a defocusing quadrupole and then a 
second non-focusing region (Focusing 0 Defocusing 0). 

The net focusing effect of an alternating gradient system may be demonstrated by calculating the 
transfer matrix for a quadrupole doublet (i.e. a defocusing quadrupole, drift space, focusing 

quadrupole arrangement) when approximating a quadrupole to a thin lens of length I, with focal 

i 

length given by Kl = -: 



M 



doublet 



1 o 

i 

7 1 



1 L 
.0 1J 



1 0 
1 

i7 1 



L 
1 



L 

1+ 7. 



Where M D = is the transfer matrix for a drift space of length, L. 



f assumes that defocusing quadrupole is equivalent to the focusing quadrupole rotated 
through 90° 



f 2 

The resulting doublet matrix represents a thin lens with focal length — > 0 18 . 



Edwards,D. & Syphers,M.(1993) An Introduction to the Physics of High Energy Accelerators 
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QF QD QF QD QF QD QF QD QF QD 
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Figure 7: A number of FODO cells, particles with non-zero divergence are confined to the 

accelerator. 

The Particle Action 

An ideal particle may be regarded as one that follows the perfect path around an accelerator (i.e. 
does not require focusing); thinking of a particle beam as consisting of only such particles is 
unrealistic. In fact, at the point of injection a beam is made up of particles with a range of 
trajectories spread around the trajectory which would lead to an ideal orbit. 

Taking into account the accelerator coordinates defined earlier, an individual particles position and 
motion in relation to that of an ideal particle may be described by its horizontal (x) and vertical (z) 

/ dx\ ( dz\ 

displacement as well as its horizontal he' = —J and vertical ( z r = — \ divergence from the ideal 

path. The focusing effect of a FODO cell upon a particle with non-zero initial amplitude or divergence 
leads to a particle exhibiting transverse oscillations around the ideal path. 

A single particle travelling through an accelerator will trace out an ellipse when the phase space 
coordinates x and x' (or z and z') are tracked through a number of FODO cells (at a specific point 
within each cell), the area of this ellipse will remain constant as the particle progresses through the 

accelerator (at a fixed value of momentum). The area of this ellipse is given by 2n] x , where J x 20 is 
the action of the particle. The relationship between the maximum transverse displacement (or 

amplitude) and action is given by Amplitude = ^2J x f3 x , (where f3 x is given definition in the 
section on tune). 



Brandt, D. (2004) Accelerator basics. Available from: 
http://cas.web.cern.ch/CAS/Warrington/PDF/Brandt5.pdf 

20 

Alexander,!.. & Fetter,J.D.W.(2003) Theoretical Mechanics of Particles and Continua 
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Adiabatic Damping 



Liouville's theorem states that the area of an ellipse in phase space, xp x (where p x is the transverse 
momentum), will be conserved regardless of magnetic beam bending or focusing as the beam 
progresses along the accelerator path. 

Liouville's theorem: Area = $>p x dx = constant 

The area, 2n] Xf expressed in terms of the phase space coordinates xx' , will shrink with acceleration. 
This observation is referred to as adiabatic damping. 

A conserved area in terms of ] x may be found by expressing the transverse momentum as: 

dx dx ds 

Px = Y™o "77 = ym 0 — — = p s x 
at as at 

Where p s is the longitudinal momentum. 

The conserved area in the xx' plane is therefore given by: 



PsYs 9 x dx = constant 



as p s « (3 sYs 

IP 

Where /3 S = with v s the longitudinal velocity and c the speed of light 



and K= \~ n 2 



1-0 



S 



The action normalised with respect to momentum is: 

Jx,norm 



18 



Chris Edmonds 



A Small Accelerator for Treating Cancer 



Tune 

The betatron tune is defined as the number of complete transverse oscillations that occur each 
accelerator turn. For an alternating gradient accelerator the betatron tune may be defined as: 

1 L ds 

Where Q is the betatron tune. 
r ds 

and <p -j^-^ is the change in phase of the transverse oscillation over a path length ds. 

The beta function, /?(s), (also referred to as the amplitude function) defines an envelope around all 
trajectories of the particles travelling through a FODO lattice. The beta function, which is determined 
by the quadrupoles, varies as function of the longitudinal position, s, and for a FODO cell tends 
towards minima at defocusing quadrupoles and maxima at focusing quadrupoles (this expresses the 
net focusing of the FODO cell). 

In the case of a constant gradient accelerator, the beta function is no longer dependent upon the 
longitudinal position, and so the tune may be given as 21 : 

1 fds_R 
Where R is the radius of the particle orbit. 

For synchrotron design it is important that the tune does not take integer or simple fraction values 
as at these points error in multipole fields (as expressed in the Taylor expansion of the magnetic 
rigidity) can drive resonance. 



Wilson, E. (2001) An introduction to particle accelerators. 
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Effects of Quadrupole Alignment Errors 

A quadrupole misaligned from the beam axis (x = z = 0) of an accelerator by (Ax,Az) (as shown in 
figure 5) will lead to particles experiencing a dipole type field error. The strength of the field is 
related to the quadrupole displacement by: 




Figure 8: An off centred quadrupole. 

For a short quadrupole (longitudinal length I) the effect of the dipole field errors influence may be 
approximated to a single angular 'kick' (with no instantaneous change in amplitude) at the midpoint 
(1/2) of the quadrupole. The nature of this kick is given by: 

Ax'\ = q/AB z \ 

Az'y p VabJ 



Wille, K. (2000) The Physics of Particle Accelerators an Introduction. 
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Figure 9: Showing a closed orbit distortion due to a dipole field error. 

Considering a particle with zero initial offset of amplitude and trajectory relative to the ideal beam 
path, (x, x') = (0,0), the amplitude and trajectory just after a kick is given by (x, x') = (0, Ax')- This 
new trajectory will, as the particle progresses longitudinally, give an increase in the transverse 
displacement, and with subsequent focusing lead to transverse oscillations. 

If the betatron tune is of an integer value, then the dipole field error kick due to the misalignment of 
a quadrupole will add coherently each accelerator turn. This means that over a number of turns the 
amplitude of the transverse oscillation will be increased, and for a sufficiently large number of turns 
beam particles will be incident upon the accelerator wall. 

However, if the betatron tune is of a non-integer value, then the dipole field error kicks of each 
accelerator turn will add incoherently and over a number of turns there will be no net change in the 
transverse oscillation amplitude. 



Wolski, A. (2006) Effects of Linear Imperfections. Linear Dynamics lecture notes. Available at 
http://pcwww.liv.ac.uk/~awolski/Teaching/LinearDvnamics/LinearDynamics-LecturelO.pdf 
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Figure 10: For integer tune (left) the amplitude grows upon each accelerator turn, for half integer 

tune (right) there is no net amplitude growth over a number of turns. 



Lee, S.Y. (2004) Accelerator physics. 



22 



Chris Edmonds 



A Small Accelerator for Treating Cancer 



Section III - Simulation 

Aims 

• To develop a computer model, based on the EMMA accelerator, which is able to simulate 
the increase in transverse amplitude due to integer tune crossing. 

o To construct the model on the premise that an alternating gradient accelerator may 
be approximated to a constant gradient accelerator for the purpose of investigating 
transverse amplitude growth. 

• To use the developed model to examine the effects of acceleration rate transverse 
amplitude growth. 

• To compare the output of this model to the output which uses a 'full physics' approach to 
describing the transverse motion. 

N.B. 'full physics' refers to a model for which the behaviour of a particles transverse motion at a 
given time is dependent upon the accelerator component in which the particle is located (i.e. the 
lattice is not approximated to a constant gradient). The exact model to be used for comparison is 
that presented by S.Machida in the paper 'Resonance Crossing and Dynamic Aperture in Non-scaling 
Fixed Field Alternating Gradient Accelerators', an outline of which may be found in the appendix. 
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The EMMA Lattice 

EMMA is the proof of principle linear NS-FFAG accelerator upon which this simulation shall be based. 




Injection 



Figure 11: Schematic of the main EMMA ring as well as injection component from ALICE and 

extraction line. 

In considering the guidance and focusing components of EMMA the accelerator ring may be divided 
into 42 periodic FODO cells. The quadrupole magnets within EMMA are combined function, with each 
magnet aligned off centre in order to provide a dipole field component for bending particles around 
the ring, and a quadrupole field component for focusing a divergent beam. 



The values used in the simulation to define the lattice may be found in the appendix. 



Conform. (2007) EMMA layout phase 2. Available at: http://www.conform.ac.uk/documents/emma/dwg- 
drawings/207-10324-C-EMMA-LAYOUT-PHASE-2-colour.jpg 
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Model Approximations 

Some of the approximations to be made in the production of the model of EMMA are: 

1. The transverse motion of the NS-FFAG may be treated approximated to that of an 
accelerator with a constant field gradient. 

2. The Transverse Motion Along the z Axis may be Modelled Independently. 

3. Electrons are accelerated from 10.75MeV to 20.75MeV, where acceleration is represented 
by a constant energy increase at each cell. 

4. The time take to complete one orbit is constant through acceleration. 

5. The lognitudinal velocity is equal to the speed of light. 



1. Constant Gradient Aprroximate Model 

The transverse motion of a particle travelling through an accelerator may be described by Hill's 
equation: 

g + K(s)z = 0 [3.1] 
and 

g + ff(s)* = 0 

Where K is dependent upon which component is being traversed within an alternating gradient 
accelerator, and is in the form of equation [2.2] for a quadrupole. 

The repeating cell arrangement of EMMA means that K is periodic, so that: 

K(s + c)= K(s) [3.2] 

Where here c may represent the length of a cell. 

A primary aim of this project is to explore the notion that for a lattice with known tune values, the 
transverse motion may be described approximately by the equation: 

g + o> 2 z = 0 [3.3] 

Where 

2nQ z 

a) = — - — 
T 

1 o 
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With Q z being the betatron tune for the vertical axis 

T 0 being theorbit time for a particle in the accelerator. 

This approximation effectively treats the focusing effect of an alternating gradient lattice as a 
constant gradient spanning the enitre circumference of the accelerator, and describes particle 
motion as sinusoidal around the ideal beam path (in the fashion of a simple harmonic oscillator), 
rather than as the more complex motion shown in figure 7. The constant gradient representation 
would see equation [3.1] interpreted as: 



i.e. K is still dependent upon the momentum of the particle, but no longer dependent upon the 
longitudinal position V. 

By rearranging these simple harmonic motion type equations [3.3] and [3.4], it can be seen that the 
restoring force may be taken as 



Relating through F = ma, gives the relationship between K and co as: 

K = moo 2 

2. The Transverse Motion Along the z Axis may be Modelled Independently 

This is the approximation that amplitude growth in the direction of the z axis may be simulated 
without consideration of the amplitude growth along the x axis. This leads to a simplified model, and 
a reduction in the amount of computer time required for simulations. This is may be a valid 
approximation if the average velocity along the x axis is a fraction of the longitidunal (5 axis) 
velocity. It is the transverse motion in the vertical (z) axis that is to be simulated as for a FFAG this 
has the smaller aperture, also the behaviour of transverse motion along the vertical axis should be 
representitive of transverse motion for the horizontal axis. 




+ Kz = 0 [3.4] 




= —Kz 



and the acceleration as 
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3. Energy Range of 10.75 MeV to 20.75 MeV with Acceleration Represented by a Constant 
Energy Increase at Each Cell 

The acceleration range of 10.75 MeV to 20.75 MeV differs from the acceleration range of EMMA 
(10.5 MeV to 20.5 MeV), but is in keeping with the 'real physics' model with which results from this 
simulation will be compared. A small increase of 0.25 MeV in injection energy keeps the initial 
horizontal and vertical tune away from integer values. 

The representation of acceleration as a constant energy increase at each cell (rather than simulating 
the effects of rf cavities every other cell) is a simplification which should result in a reduction in the 
amount of time each simulation takes to run. Once again, this is in keeping with the 'real physics' 
model, for which, this representation of acceleration also eliminates the effects due to longitudinal 
and transverse coupling 26 . 

4. Constant Time of Flight 

For the purpose of this simulation, the time taken for an electron to complete one turn in EMMA will 
be regarded as being constant with energy, and is to be defined by the circumference of EMMA 
divided by the speed of light (i.e. independent of the amplitude of transverse motion). 
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Graph 3: Time of flight is approximately isochronous within the range of acceleration. 



26 For a model in which the time taken for one accelerator turn is dependent upon the transverse 
amplitude of a particle, it is possible that a particle with a sufficiently large transverse amplitude will 
enter into the decelerating phase of the rf. Representing acceleration as a constant increase in 
energy per cell isolates this problem from that of particle loss due to resonance. 



Berg, J. (2007) The EMMA Experiment. 
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The NS-FFAG accelerator lattice is designed so as to give a parabolic variation in orbit length during 
acceleration (with the minimum of the parabola located at the median particle energy). Referring to 
graph 5, the range in time of flight deviation per cell, 3ps, is small compared to the time of flight per 
cell, which is approximately 1. 3 ns (for a particle with orbit length of 16.57 m and velocity close to 
the speed of light). 

When considering the time of flight dependency on transverse amplitude, the relative differences 
between the wavelength of the transverse oscillations and the transverse amplitude may be taken 
into account. The shortest wavelength (A) of vertical transverse oscillation within EMMA is 1.3m; the 
maximum vertical transverse amplitude (A) (limited by the vertical aperture in the focusing 
quadrupole) is 8.9mm. An estimation of the maximum difference in path length between a zero 
amplitude particle and a high amplitude particle is obtained by approximation of the path length 
between the maximum and minimum of the oscillation to being the hypotenuse of a triangle with 
sides 2A and A/2: 




Figure 12: Estimation of path length for a high transverse amplitude particle in EMMA. 

This method shows an increase in path length of 0.04% of the straight through path for a high 
amplitude particle. 
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5. Longitudinal Velocity Equal to the Speed of Light 

For the electron energy range of 10.75 MeV to 20.75 MeV there is a 8 x 10~ 4 c (where c is the 
speed of light) increase in speed, and a maximum deviation in speed from the speed of light of 
1 x 10" 3 c (a difference of 0.1%). 

0.9998 




E (MeV) 



Graph 4: Beta Vs Energy for EMMA's acceleration range 

The speed of the particle is given by <sjv x 2 + v z 2 + v s 2 , the assumption is made that 

v s » V^x 2 + v z 2 (from this it follows that the Lorentz force due to transverse velocity will be 
relatively small). 
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A solution to the simple harmonic motion equation of [3.3] comes in the form of: 

z = Acos(a)t + cp) 

Which may be rewritten using the trigonometric identity 
cos(x + y) = cos(x)cos(y) — sin(x)sin(y) to give: 
z = Acos(a)t) + Bsin(a)t) [3.5] 
The velocity of a particle at a time, t, is then given by: 

dz 

v — — — —(x)Asin{(jdt) + a)Bcos(a)t) [3.6] 



and the coefficients A and B are found by setting t = 0. 
At t = 0, z t = A 

(where the subscript T denotes initial) 
Substituting into [3.5] and [3.6] gives 

z = ZjCos(a)t) + — sin(a)t) 



v = —ojzisin^ajt) + ViCos(a)t) 

Which may be used to create a transfer matrix for transverse displacement and velocity through a 
time t: 



This transfer matrix will be used to estimate (under the constant gradient approximation) the 
transverse amplitude and velocity of a particle at intervals consistent with the location of individual 
quadrupoles. 

The transverse amplitude is then given by: 



All particles within the simulation are started with 0mm initial amplitude. 



OJ 




cos(a)t) 
—ti)sin(a)t) 




[3.7] 
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The main features of the Matlab code for doing this are annotated below (in red): 

for cl = 1 : (CT-1) , 

1. Loop runs over the total number of cells traversed during acceleration, i.e. for 
acceleration over 100 turns, CT=100*42 

energy = Estart + cl*deltaE; 

tune = (0 . 0549*energy A 2) - (2 . 4572*energy ) + 32.4479; 
freq = tune / TO; 
wb = 2*pi*freq; 

2. Energy and tune (frequency) parameters unique to each cell. 

dvF = (F(c,3) * Light Speed * Q (2, 1) ) / (energy*le6*qe) ; %transverse 
velocity kick at focussing quad 

dvD = (F (c, 4) * LightSpeed * Q (2, 2) ) / (energy*le6*qe) ; ^transverse 
velocity kick at defocussing quad 

3. Gives transverse focusing and defocusing velocity kicks based upon misalignment of 
current cell. 

vectorl = [xl(cl) 

vl (cl) ] ; 

4. Initial transverse displacement and velocity, based upon displacement and velocity at 
midpoint of defocusing quadrupole in previous cell (after kick), or if cl=l, upon the 
injection conditions. 

M = Ml (wb, TDF) ; 

5. Gives transfer matrix specific to the angular frequency of the particle and to the time 
taken for the particle to travel from the midpoint of the defocusing quadrupole in the 
previous cell, to the midpoint of the focusing quadrupole in the current cell. 

vector2 = M * vectorl; %travel from middle previous defocussing quad to 
middle of following focussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + dvF; %adds kick from horizontally focussing quad 
to vector 

6. Estimates transverse displacement and velocity (based upon constant gradient 
approximation) of particle at midpoint of focusing quadrupole, and adds transverse 
velocity kick. 

M = Ml (wb, TFD) ; 

vector2 = M * vector2; %travel from previous focussing quad to following 
defocussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + dvD; %adds kick from horizontally defocussing 
quad to vector 

7. Carries out steps 5 and 6 for the defocusing quadrupole (based on time taken for 
particle to travel from midpoint of focusing quadrupole to midpoint of defocusing 
quadrupole. 

xl(cl+l) = vector2(l); 
vl(cl+l) = vector2(2); 

8. Sets initial transverse displacement and velocity for next cell. 

Al(cl+1) = sqrt (xl (cl+1) A 2 + (vl(cl+l) A 2 / wb A 2)); 
betal(cl+l) = L/ (2*pi*tune) ; 

9. Calculates amplitude and betatron function. 

c = c+1; 
if c>42 
c=l; 
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end 

10. For loop counts through cells, and repeats after 42 (i.e. one revolution), 



end 



Plots of single particle amplitude Vs tune are presented in the form of ^2f z Vs tune to aid 
comparison with the real physics model. 

The transverse displacement, z, is related to ^2T Z by: 
z = ^2] z (3 z coscp z 

z is a maximum when coscp z = 1, which gives the amplitude of the oscillation along the z axis as: 

A z — V ^Jzfiz 
A, 



With the approximation of an alternating gradient accelerator to a constant gradient accelerator the 
beta function, /? z , is longitudinally invariant. ji z , given as a function of particle energy is: 

R 



where R is the radius of the particle orbit 
Q Z (,E) is the vertical tune at energy 'E'. 



The plot of y]2J z Vs tune for no perturbing influence and a starting displacement of 1mm is as 
follows, the effects due to adiabatic damping are not applicable (which is contrary to the output of 
the real physics model): 
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Graph 5: Start amplitude of 1mm and standard deviation of alignment error magnitudes of 0m. 
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Introduction of Quadrupole Alignment Errors 

The alignment errors to be considered are for the offsetting of the cell position (the focusing and 
defocusing quadrupole within cell are equally misaligned) into the x-z plane, the magnitude of the 
errors are assumed to be of a Gaussian nature with a 2a cut off. 

The Gaussian probability density function is given by: 

1 /-(AR-|i) 2 \ 
P (AR)= — =exp \ 2 W 

Where P(AR) is the probability of an error of a value AR 
\i is the mean 
a 2 is the variance 

The magnitude of the quadrupole alignment error is produced using the Matlab 'randn' function, 
values obtained falling outside of the 2o range are reevaluated: 
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Graph 6: Blue line showing the Gaussian probability distribution for [i=0 and o=2xl0 7 m. Red 
histogram showing the result of lxlO 6 trials based on the same parameters as for the probability 

distribution (blue line), and with outcomes of \AR\ > 2a re-trialled. 
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A direction is then assigned to the error using the Matlab 'rand' function (which generates a pseudo- 
random number in the interval 0-1), through 6 = rand * 2n. 

AR X = AR cosd 
AR Z = AR sind 

Figure 13: Alignment is offset 
at a random angle into the x-z 




The distribution of quadrupole alignment errors along the z axis (the axis for which amplitude 

growth will be considered) is as follows: 
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Graph 7: Showing the distribution of errors in quadrupole alignment along the z axis. Generated by 
lxlO 6 trials. 
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When studying the amplitude growth of a single particle a 42 element array is produced to contain 
the magnitudes of cell alignment errors. Random number seeds (to determine the set of pseudo- 
random numbers utilised) are used in order to provide direct comparison between different rates of 

acceleration. Below is the set of cell alignments used to produce graph 10 (^2f z Ms tune): 

x 10" 5 

6 i , , , , , , , , , 




Graph 8: Quadrupole alignment errors in z direction as utilised in simulation for graphs Seed 

number 788. 



The effect of a misaligned quadrupole may be seen as a kick in the transverse velocity, which is given 
by the Lorentz force: 

F = q(v x B) 

For an error in the quadrupole alignment along the z axis, the Lorentz force is given by: 
F z = qv s ^Az (3.8) 

Which leads to the acceleration of the particle along the z axis: 
F z = m^ (3.9) 

z dt v ' 

E 

Where m is the relativistic mass, m = ym 0 = — 

C 

The transverse velocity kick is defined by equating (3.8) and (3.9): 



35 



Chris Edmonds 



A Small Accelerator for Treating Cancer 



F z QV S dB x 

dv 7 — — dt = —Azdt 

m m dz 

And is applied at the midpoint of a misaligned quadrupole: 

v z j = v z ,i + dv z 
Where 

v z j is the transverse velocity after the kick 
v z i is the transverse velocity before the kick 
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Effect of the Standard Deviation of Quadrupole Misalignment on 
Amplitude Growth 



It can be expected that the extent of transverse amplitude growth at integer tune values is related to 
the magnitude of the standard deviation of quadrupole alignment errors. In order to test this 
proposition, simulations were carried out for quadrupole alignment errors of standard deviations 
ranging from 0.2[im to lOO^m (at intervals of 0.2\im). For each value of standard deviation a unique 
pattern of quadrupole alignment errors was created, a particle (with zero initial transverse 
amplitude) was then accelerated from injection energy to extraction energy and the maximum 
transverse amplitude throughout acceleration recorded. A total five hundred error patterns was 
selected so as to reduce statistical errors when finding a line of best fit. For all simulations (unless 
otherwise stated), particles start with zero transverse displacement and velocity, any non-zero 
amplitude found during acceleration is therefore as a result of growth due to quadrupole 
misalignments. 

A plot of maximum amplitude Vs standard deviation was found to be as follows: 
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Graph 9: Plot of maximum amplitude Vs standard deviation for acceleration over 100 turns, line of 
best fit generated using the Matlab polyfit function. 
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Rate of Acceleration 



The number of accelerator turns for which resonant conditions may lead to amplitude growth is 
dependent on the number of turns over which accelertion occurs, with slow acceleration rates 
leading to stronger amplitude growth at integer tune values. Insight into the importance of the rate 



of acceleration may be gained visually bycomparing plots of y 2] z Vs Tune for regimes identical in all 
but the number of turns over which acceleration occurs: 
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Graph 10: Top acceleration through 100 turns, bottom acceleration through 1000 turns. Both plots 
completed using the same pattern of quadrupole alignment errors (standard deviation of 30|im). 

In order to gain a quantitative appreciation of the effects of the acceleration rate, the simulations 
completed in order to produce graph 9 (maximum amplitude Vs standard deviation) are repeated for 
a range of acceleration rates (acceleration over 100, 200, 300, 400, 500, 750, 1000, 1500 and 2000 
turns). Below is a plot showing maximum amplitude Vs standard deviation for samples of 
acceleration over 100 and 1000 turns: 
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Graph 11: Transverse amplitude Vs Standard deviation of alignment error for acceleration through 
100 turns and 1000 turns, completed for 500 different values of standard deviation (each 

representing a different pattern of alignments). 



For each rate of acceleration a line of best fit (in the form of y = mx + c) for the data of maximum 
amplitude Vs standard deviation is produced using the Matlab polyfit function. A graph of the 



gradient of each line of best fit, 



max 



^ , is then plotted as a function of the rate of 



acceleration. 



The error on each value of 



max 



^ a is based upon the distribution of the points of j2J z 



max 



around the line of best fit, and is calculated by 28 : 



®m 




(iV-2) 



N 



Where, 



N is the number of points 



Press, W.H. et al. (1995) Numerical Recipes in C 
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The x 2 value is given by: 



X 



2 =i 



N x 2 

1 hi - (rrtxt + c) 



1=1 



Where, 

yi is an individual measured point (in this case a value for maximum amplitude obtained 
through simulation) 

mxi + c is the value of y t as determined by the line of best fit 

&i is the error on the individual point 

In determining the error on the gradient each value for maximum amplitude is treated with equal 
weighting (<T/ = 1). 



max 



^ plotted as a function of the rate of acceleration: 
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Graph 12: Amplitude growth as a function of the rate of acceleration. 
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is the inverse of square root of the tune change per turn, and is chosen so that data 



dT 



points sit on a straight line (in keeping with figure 4 of Machida (2008)). In calculating ^/^j the 



tune is described as varying linearly with the number of turns, so that: 



dQ I Maximum Tune -Minimum Tune 

'dT ~ Number of Turns 



The gradient describing a line of best fit to the points of graph 12 is found as that which minimises 
the x 2 function I dm ~ ® where m is the gradient^ 29 . 



m = 



SS X y S X Sy 



S 
A 



Where the individual terms are defined as: 



c _ V N 1 



*xy — Zji=i : 



c = y N c — v^v y±_ 



^— ^xx S x S xx — Yii=i 



The gradient obtained for the amplitude growth as a function of acceleration rate was: 



A l /Act 



= (38 + 1 )m _1 /2 (3% error) 



\dQj 



dT 



With the fit giving a x 2 P er number of degrees of freedom as 1.27. 
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Press, W.H. et al. (1995) Numerical Recipes in C 
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The Physical Aperture 



The physical aperture was defined as the minimum initial amplitude (for a particular rate of 
acceleration and distribution of quadrupole alignment errors) for which a particle was not 
accelerated to the extraction energy as a result of the particle being incident upon an accelerator 
component. 

Within the simulation particles with ^2f z of more than 19.3rara~ '2 were considered to have been 
lost, this is consistent with the transverse amplitude of a particle being more than half the aperture 
of the smallest aperture component (this is for the focusing quadrupole which has a vertical 
aperture of 17.8mm). 

Five hundred different values of the standard deviation of quadrupole misalignment are evaluated 
for two different rates of acceleration (through 100 and 1000 turns). 

For each combination of standard deviation and acceleration rate, the simulation runs trials with the 
start amplitude of particles increasing from 0mm at intervals of 0.1mm until the criteria of 

J2J Z > 19.3mm~ '2 is met at some point during an acceleration cycle. 
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Graph 13: Scatter plot of physical aperture Vs Standard deviation for acceleration through 100 and 

1000 turns. 

The physical aperture is displayed in units of (u ITI rad). 

To give greater clarity when visually interpreting the results displayed in graph 13, averages of the 
physical aperture over intervals of 15 standard deviations are found and plotted: 
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Figure 14: Physical Aperture Vs Standard Deviation (physical aperture averaged at intervals of 15 

standard deviation values). 
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Section IV - Analysis and Conclusion 
Analysis and Comparison With a Real Physics Model 



Where possible, the parameters of the simulation were set so as to match those used in the model 
presented by Machida (2008) in the paper 'Resonance Crossing and Dynamic Aperture in Non-scaling 
Fixed Field Alternating Gradient', this was to allow for the evaluation of the approximation that an 
alternating gradient accelerator may be represented as having a constant gradient when 
investigating transverse motion behaviour. 

Single Particle ^f2U Vs Tune 

Plots were produced for the amplitude growth of a single particle as a function of tune (and so 
energy). The graphs below represent different distributions and standard deviations of alignment 
errors, and so the meaning of direct comparison is limited. Both plots have been smoothed using 
Savitzky-Golay filtering (see appendix): 




Graph 14: Amplitude Vs tune for a real physics model (top) and the constant gradient approximation 
model (bottom), both with acceleration through 1000 turns. 



Figure 2 of 'Resonance Crossing and Dynamic Aperture in Non-scaling Fixed Field Alternating Gradient Accelerators' The 
axis referred to as 'y' in this plot is equivalent to the axis referred to as 'z' throughout this document. 
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The amplitude values for the 'real physics' model have been normalised with respect to /?y for this 
plot of amplitude growth. Also, for the 'real physics' model particles were injected at an amplitude 
which was consistent with the displacement of the phase space origin for a lattice with errors, from 
the phase space origin without errors, whereas for the constant amplitude model particles were 
injected with zero amplitude. This finite initial amplitude for the 'real physics' model is assumed to 
be small when compared to the amplitude growth driven by integer tune crossing. 

Features consistent with dipole field error driven resonance at integer tune crossing may be seen for 
both models, these are the sharp changes in amplitude in the vicinity of integer tune values and a 
relatively small variation in amplitude when the tune is away from integer values. 
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Effect of the Standard Deviation of Quadrupole Misalignment on Amplitude Growth 



By plotting graphs of y 2] z Vs Tune for conditions which may be differentiated only by the standard 
deviation of the alignment error (i.e. same rate of acceleration, same transverse start amplitude and 
phase and same error patterns 31 ) the following may be obtained: 

x 10" 3 




1 1 1 1 _L_ L 

12 11 10 9 8 7 6 

Tune 



Graph 15: Top quadrupole misalignment distribution with standard deviation of 10 |im, bottom 
with standard deviation of 50|am. Error patterns identical (i.e. same seed number used for errors). 

In essence the plot for errors with standard deviation 50|im is an amplification of the plot obtained 
for the lO^im distribution. This seems unrealistic, and is potentially as a result of defining the 
longitudinal velocity as being constant (particularly in the sense of being independent of transverse 
velocity). In order to gain better understanding of the implications of this approximation, it would be 
advantageous to develop the model further and give a better representation of a changing 
longitudinal velocity. 



3 1 

The code for generating quadrupole alignment errors is of the form 'error= u.+ a*randn(l)'. Altering the 
standard deviation of the guassian which describes the distribution of quadrupole alignment error essentially 
scales the errors accordingly. 
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Effect of the Rate of Acceleration 



When making comparisons between the 'real physics' and constant gradient approximated model 



for plots of the maximum amplitude through acceleration (yj2J z>max ), it is worth noting that 

Machida seems to suggest that the values of ^2f z for the 'real physics' plot have not been 
normalised with respect to /?y (this is contrary to the 'real physics' plot presented in graph 14). The 



effect of adiabatic damping in this situation will be to give a reduction on the values of ^2] z 
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Graph 16: Maximum amplitude Vs standard deviation of quadrupole alignment errors for 'real 
physics' model (top) and constant gradient approximated model (bottom). 

Visual inspection of the plots presented in graph 16 suggests that the rate with which maximum 
amplitude is found to increase as a function of standard deviation is similar for both models. 



Figure 3 of 'Resonance Crossing and Dynamic Aperture in Non-scaling Fixed Field Alternating Gradient Accelerators' 
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Comparing graphs of the rate of the maximum amplitude as a function of the standard deviation of 
quadrupole alignment errors Vs the rate of acceleration confirms this: 
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Graph 17: Amplitude growth as a function of acceleration rate for 'real physics' (top) and constant 

gradient approximation (bottom). 
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The gradients of the line of best fit for graph 17 are: 



\ 



,max / 
A l /Act 



i 



- (50 + 10 )m ^2 (20% error) for the 'real physics' model 



IdQj 



dT 



\ 



max , 



Act 



/ 



= (38 + 1 )ra ^2 (3% error) for the constant gradient approximate model 



\&Qj 



Using a standard consistency check (\x ± — x 2 \ < + a 2 2 ) for two measurements (x + cr) of 

the same quantity the values are found to be consistent: 



50 - 381 <3 /10 2 + 1 2 2 



However, it may be questionable as to whether the gradients presented are measurements of the 
same quantity. This is dependent upon the interpretation of the effects of adiabatic damping on the 
gradient obtained from the plot for the 'real physics' model being correct. Further analysis would be 
needed in order to confirm the consistency suggested above. 
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Physical Aperture 

Transverse amplitude growth was considered in terms of the apertures of the accelerator, and the 
minimum initial start amplitude for which a particle would be lost was found for a range of 
quadrupole error distributions: 




0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 



Graph 18: Physical aperture Vs standard deviation of alignment errors as given by the constant 

gradient approximate model. 

For the slower rate of acceleration (through 1000 turns) the initial amplitude for which a particle 
may be accelerated until the extraction energy falls quickly with increasing quadrupole alignment 
errors. 

For the 'real physics' model a plot of dynamic aperture Vs the standard deviation of alignment error 
reveals similar behaviour: 
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1 a of alignment error [fim] 3 4 

Graph 19: Dynamic aperture Vs Standard deviation of quadrupole alignment errors as found with 

the 'real physics' model. 

Machida defines the dynamic aperture as being the 'maximum amplitude of survival particles', 
where a survival particle is one for which the transverse amplitude does not grow to beyond the 
stable region in transverse phase space. As with the simulation of physical aperture (for the constant 
gradient approximation), the above plot is obtained by locating the maximum 35 initial amplitude for 
which a certain value of amplitude is not met throughout acceleration. 

What is not immediately clear, is the reason behind the difference in dynamic aperture for 
acceleration through 100 turns and 1000 turns when there are no quadrupole alignment errors 
present (for the physical aperture plot the lines for both rates of acceleration converge with no 
alignment errors). 



Figure 5 of 'Resonance Crossing and Dynamic Aperture in Non-scaling Fixed Field Alternating Gradient Accelerators 7 
35 For the relative size of interval the maximum initial amplitude for which a condition is not met may be 
approximated to the minimum initial amplitude for which a condition is met. 
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Conclusion 

It is difficult to give an evaluation of the quantitive agreement between models without 
understanding the exact parameters of both; however these initial findings suggest that there is 
consistency. 

The approximation of alternating gradient focusing to constant gradient focusing was found to be 
effective in demonstrating the fundamental behaviour of transverse amplitude growth when 
resonance at integer tune values is driven by dipole field errors; With the ability of a linear NS-FFAG 
to accelerate a particle to an extraction energy being dependent upon the rate of acceleration and 
the magnitude of the quadrupole alignment errors within the accelerator. 

Looking at figure 7, it can be seen that the relationship between the phase of the transverse motion 
and the location of quadrupoles is dependent upon the properties of the lattice. This is something 
that has not been accounted for a within the constant gradient approximated model, and the effects 
of which, I feel, are worthy of further investigation. 

Neither model explores the effect of resonances as a result of tune variation within EMMA fully, with 
only dipole field errors considered. Higher order errors, such as quadrupole field errors, will 
introduce resonances at non integer tune values (when the tune value is a simple fraction, e.g. half 
integer) 36 . A more complete model could investigate these effects of these too. 

In terms of the suitability of the NS-FFAG for proton therapy applications, the attraction to a 
relatively low cost combination of the intensity of a cyclotron and the flexibility of a synchrotron is 
clear. For the linear NS-FFAG studied, the nature of amplitude growth with relation to the standard 
deviation of quadrupole alignment errors and rate acceleration may mean that a non-linear NS-FFAG 
lattice design 37 must be considered for proton acceleration(Machida(2008)). 



36 

Brandt,D (2004) Accelerator Basics presentation. Available at 
http://cas.web.cern.ch/CAS/Warrington/PDF/Brandt5.pdf 

37 Which may limit the tune excursion in order to exclude integer tune values. 
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Section V - Appendix 

Summary of the article: 

Resonance Crossing and Dynamic Aperture in 
Non-scaling Fixed Field Alternating Gradient Accelerators 



The article, written by Shinji Machida, describes an investigation into the effect of quadruple 
alignment errors and the rate at which particle acceleration occurs upon the amplitude of a particle 
for a linear nsFFAG. Within the abstract for the article Machida reports that for a typical error of 
lO^im rms, that there is 'practically no dynamic aperture' for slow acceleration (through 1000 turns). 

The following offers a brief summarisation of the article, with the aim of presenting Machida's 
simulation results (which shall act as a comparison for the output of the simulation being developed 
as part of the Phys498 project). 

JjJntrP.ductip.n 

Machida begins by introducing the concept of accelerating muons (to energies of 20 or 50 GeV) 
through the use of a scaling (s) FFAG. The betatron tune for a sFFAG remains constant for both 
planes through acceleration, however, for a nsFFAG this is not the case as tune varies with inverse 
proportionality to energy (see graph below). As a particle is accelerated the tune crosses through 
integer values, and it is in these regions that kicks from misaligned magnets lead to dramatic growth 
in particle amplitude (and the eventual loss of the particle to the accelerator wall). 




12 14 16 IS 20 
momentum [MeVte] 

Graph 20: Tune Vs Momentum 

The muon lifetime (2.2\xs) means that rapid acceleration must be utilised in order to reach the 
energies quoted above (with acceleration occurring over 10 to 20 turns), in this instance there is no 
destructive build up in amplitude as a result of traversing integer tunes. For other anticipated 
applications of nsFFAG's, such as for particle therapy, priority shifts from rapid acceleration (proton 



Machida, S. (2008) Resonance Crossing and Dynamic Aperture in Non-scaling Fixed Field Alternating Gradient 
Accelerators 
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lifetime many orders of magnitude higher than that of the muon) to reducing costs (accelerators 
with lower rf voltages). 

It is for these slower rates of acceleration that resonance crossing becomes problematic, and, for 
this reason, Machida raises the importance of determining the feasibility of such acceleration 
regimes. 

J J .■ . Sim. u I [a t i.o n Se t up, 
A. Acce [era to_r_M pd e I 

The EMMA lattice is used as a basis for simulation. 
Several assumptions were made: 

• Particle energy increased by a constant amount at each cell (to represent acceleration). 

/. Discounts the effects of longitudinal and transverse coupling: 

For the EMMA lattice, the time take for one revolution for a particle is dependent 
upon the transverse amplitude (as larger amplitude particles follow a larger path). 
This difference in revolution time leads to energy spread in a bunch, and for a 
particle with amplitude beyond a critical magnitude, entering the deceleration phase 
of the rf will lead to the particles loss. 

ii. Leads to a reduction in processing time for the simulation. 

Hi. Acceleration occurs at each cell rather than at periods consistent with the location of 
the 19 rf cavities so as not to create additional symmetries within the lattice. 

• Simulation injection and extraction energies of 10.75 and 20.75 MeV (rather than EMMA 
values of 10.50 and 20.50 MeV). 

Acceleration begins with horizontal and vertical tune values away from integer 
quantities and results in lower initial orbit distortion. 

• Transverse amplitude growth is stimulated by quadrupole misalignment. The focussing and 
defocusing quadrupoles of a cell, which are located on the same table, are assumed to be 
misaligned in the same direction and by the same magnitude. Rotational and longitudinal 
misalignments of quadrupoles are neglected. 

_B. Partjcl e _Am p li t_u_d_e 

The value of ] y is the action for the particle in yy' plane, values will exhibit adiabatic damping. 
nv y 2 + {a y + p y y') 2 
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J J J .■ .Sj m uj a t ion. Res u Its 

A. Single P a rti c I e _B_e_h ay iour 

Machida presents the development of normalized (with respect to relativistic (3y) vertical amplitude 
with decreasing tune (i.e. increasing particle energy). Vertical motion is selected due to its similarity 
to horizontal motion and the smaller vertical aperture. 
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Graph 21: Normalised amplitude Vs tune for (a) 100 turns and (b) 1000 turns. Magnitude of 
misalignment in (a) higher than that of (b) hence larger amplitudes for (a) 

Graphs 2. a and b. demonstrate how amplitude growth is related to the rate of acceleration; with the 
slower rate of acceleration (through 1000 turns) showing dramatic amplitude changes at integer 
tune values. 
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The rate at which the maximum transverse amplitude increases with the standard deviation of the 
quadrupole misalignment error is shown to be dependent upon the rate of acceleration: 




1 G of alignment error [|im] 

Graph 22: Maximum particle amplitude (through acceleration) Vs the standard deviation of the 
alignment errors for acceleration through 100 and 1000 turns. 



The growth of maximum amplitude with respect to acceleration rate: 
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Graph 23: (Figure 4 form Machida's paper) Amplitude growth as a function of acceleration rate. 



According to Machida, graph 20 has a gradient of 50m" 1/2 5 based upon the same method of error 
calculation as utilised for the constant gradient approximate model the value with uncertainty can 
be taken as (50±10)m" 1/2 . 
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_B. D_y na rn ic Ajperture 

Machida explains how the growth of amplitude may result in a particle not being accelerated to the 
final energy and defines the dynamic aperture as being the maximum initial amplitude of 'survival' 
particles. Here a survival particle is classed as one which is not accelerated to outside the region of 
stable phase space. 




0 20 40 60 80 100 
1 a of alignment error [|im] 

Graph 24: (Figure 5 from Machida's paper) Dynamic aperture Vs standard deviation of alignment 
error, with solid lines showing average of nearby points. 

Machida notes that there is no dynamic aperture when quadrupoles have alignment errors of more 
than 15[im when acceleration occurs over 1000 turns. 
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Testing of Basic Model 

Zero amplitude does not grow without perturbation: 
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Graph 25:Start amplitude Omm and standard deviation of alignment error magnitudes of Om. 

Amplitude (mm) increases as a function of energy regardless without the effects of dipole field 
errors. Model equivalent to simple harmonic motion for a spring with constantly weakening spring 
constant k: 
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Graph 26: Start amplitude of 1mm and standard deviation of alignment error magnitudes of Om. 

Amplitude in terms of ^2] is constant with increases energy. 

The effect of a transverse velocity kick upon an oscillator can be to increase the transverse 



amplitude: 
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Graph 27: The effect of a periodic perturbing influence (at intervals of 65ns) on a 960 Mrads 
oscillator. 



Presentation of Results for ^]2] 7 Ms Tune 

Savitzky-Golay filtering is a method that smoothes noisy data whilst preserving features such as 
minima and maxima. The filter functions by fitting a polynomial of the n th degree at intervals of x 
(where n and x are variable) to data. 
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Graph 28: Top-Without Savitzky-Golay filter. Bottom-With Savitzky-Golay Filter 
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Figure 15: EMMA cell layout. 

F0D0 properties: 



Number of cells in EMMA 


42 


Long drift length 


0.210m 


Short drift length 


0.055m 


Defocusing quadrupole Length 


0.071m 


Defocusing quadrupole gradient 


-4.892Tm 1 


Focussing quadrupole length 


0.058m 


Focussing quadrupole gradient 


6.650Tm 1 


Total length 


16.573m 



Table 2: Some of the EMMA geometrical information relevant to this project. 

Quadrupole Aperture data: 





Defocusing 


Focusing 


Separation 


50mm 


Yoke lengths 


65mm 


55mm 


Inscribed radii 


51mm 


36mm 


Aperture (h) 


26.2mm 


42.3mm 


Aperture (v) 


23.4mm 


17.8mm 


Max. field (T) 


0.44 


0.48 



Table 3:focusing and defocusing quadrupole parameters. 



Conform (2007) EMMA cell layout. Available at 
http://www.astec.ac.uk/emmafiles/magnets/EMMAcelllavout.ppt 

40 Conform (2007) EMMA geometry. Available at http://www.astec.ac.uk/emmafiles/magnets/EMMA 
Geometrv.xls 
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Betatron tune: 
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Given that the betatron tune is a function of the lattice, the relationship between betatron tune and 
energy was defined in the simulation by a 2 nd order polynomial fit to the data presented by 
S.Machida. 

Q z = 0.0549£ 2 - 3.4572£ + 32.4479 
Q z = 0.0621£ 2 - 2.5827£ + 33.7552 



Where E is the energy. 
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Graph 29: Plot of vertical (blue) and horizontal (green) tune as a function of energy. As given by 

polynomial fit. 



Edgecock,R. (2007) EMMA - THE WORLD'S FIRST NON-SCALING FFAG 
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Matlab Code 

1. Matlab program (ex3multikickcurvytune.m) for producing results in the form: 



x 10 

8 r 




0 [ 1 [ ' ' L 

12 11 10 9 8 7 6 

Tune 



Main Program: 

clear all 

% to simulate the effect of 84 kicks per revolution 

% on an electron in EMMA 

% assuming uniform focusing strength 

o, _ 
o 

^Constants throughout simulation 

global Q qe LightSpeed 

NO = 1000; % no. of turns 

CT = 42*N0; %number of cells traversed 

L = 16.57; % (m) Circumference 

%Quadrupole details ( Column 1 is focussing and 2 is defocussing, 
%row 1 is field gradient and 2 is quadrupole length) 
Q = [6.7, -4.9 

58 . 782e-3, 75 . 699e-3] ; 

qe = 1.602e-19; ^electron charge 

LightSpeed = 2.998e8; % (m/s) speed of light 
TO = L/LightSpeed; %time taken for one revolution 

TDF = 0 . 274571/LightSpeed; %time taken to travel from defocussing to 
focussing quad centre 

TFD = 0 . 120023/LightSpeed; %time taken to travel from focussing to 
defocussing quad centre 

o, _ 
o 

% set matrices of before and after position and velocity 
Al = zeros (1, CT) ; 
xl = zeros (1, CT) ; % 

vl = zeros (1, CT) ; % xl (m) displacement & vl (m/s) velocity after one cell 
betal = zeros (1, CT) ; 
tunegraph = zeros (1,CT); 

o, _ 
o 

^Acceleration parameters 
%Energy range 

Estart = 10.75; %Injection energy MeV 
Estop = 20.75; ^Extraction energy MeV 

deltaE = (Estop - Estart) /CT; %energy increase per cell 
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o, _ 
o 

%dist ribut ion of misalignments and associated forces 
mu=0; %mean 

sigma=20e-6; %standard deviation 

seed = 774; %seed for repeating results with pseudo random number 
F = gauerr (mu, sigma, seed); 

o, _ 
o 

% Particle Start parameters 

tunemin = ( 0 . 054 9*Estop A 2 ) - (2 . 4572*Estop) + 32.4479; 
tunemax = ( 0 . 054 9*Estart A 2 ) - (2 . 4572*Estart ) + 32.4479; 
betal(l) = L/ ( 2 *pi*tunemax) ; % beta start value 
freq = tunemax / TO; 
wb = 2*pi*freq; 

o, _ 
o 

Al(l)=0e-3; % (m) start amplitude 

xl(l) = 1*A1(1); % (m) start displacement 

vl(l) = sqrt (wb A 2* (Al (1) A 2-xl (1) A 2) ) ; % (m/s) start velocity 

o, _ 
o 

%count ing 

c=l; %cell number 

for cl = 1 : (CT-1) , 

energy = Estart + cl*deltaE; %accelerart ion 

tune = (0 . 0549*energy A 2) - (2 . 4572*energy ) + 32.4479; 

tunegraph (cl) = tune; 

freq = tune / TO; 

wb = 2*pi*freq; 

dvF = (F (c, 3) * LightSpeed * Q (2, 1) ) / (energy*le6*qe) ; ^transverse 
velocity kick at focussing quad 

dvD = (F(c,4) * LightSpeed * Q ( 2 , 2 ) ) / (energy* le6*qe) ; %transverse 
velocity kick at defocussing quad 

vectorl = [xl(cl) 

vl (cl) ] ; 

M = Ml (wb, TDF) ; 

vector2 = M * vectorl; %travel from previous defocussing quad to 
subsequent focussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + dvF; %adds focussing quad kick to vector 

M = Ml (wb, TFD) ; 

vector2 = M * vector2; %travel from previous focussing quad to 
subsequent defocussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + dvD; %adds defocussing quad kick to vector 

xl(cl+l) = vector2(l); 
vl(cl+l) = vector2(2); 

Al(cl + 1) = sqrt (xl (cl + 1) A 2 + (vl(cl + l) A 2 / wb A 2)); 

betal(cl + l) = L/ (2*pi*tune) ; 

Al (cl + 1) = Al (cl + 1) /sqrt (betal (cl + 1) ) ; 

c = c+1; 

if c>42 %42 cells in EMMA 
c=l; 

end 
end 

o, _ 
o 
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turn = 0: (N0-1); % turn number for graph 

SA1 = sgolayf ilt (Al, 1, 301) ; %Savit zky-Golay filter 

subplot ( 2 , 1 , 1 ) 

plot (tunegraph, SA1) %plots processed data 

set (gca, ' XDir 1 , ' reverse ' , ' XLim ' , [ tunemin tunemax] ) 
xlabel ( 1 Tune ' ) 

ylabeK'SQRT (2J_{z}) [m A {l/2}]') 
subplot (2,1,2) 

plot (tunegraph, Al) %plots raw data 

xlabel ( 1 Tune ' ) 

ylabeK'SQRT (2J_z) [m A {l/2}] f ) 



Random Cell Alignment Generator (gauerr.m): 

%function to produce an array of errors 
function F = gauerr (mu, sigma, seed) 
global Q qe LightSpeed 

err = zeros (42,1); %Quad error matrix 

F = zeros (42, 4); % Quad error force matrix 

randn ( ' seed' , seed) ; 

rand ( ' seed' , seed+100) ; 

for x=l : 42 

z = 2* sigma +1; 

while z >(2*sigma) %for as long as the error magnitude is larger than 
2*sigma cut off 

errmag=mu+ sigma . * randn ( 1 ) ; 

z = sqrt (err (x) A 2) ; 

end 

theta = rand * 2*pi; 

err (x, 1) = errmag * cos (theta) ; 

err (x, 2) = errmag * sin (theta); % gives x and y components of 
alignment error 

F (x, 1 ) =qe*Light Speed*Q ( 1 , 1 ) *err(x,l); %horizontal focussing force 

F (x, 2 ) =qe*Light Speed*Q ( 1 , 2 ) *err(x,l); % horizontal defocussing force 

F (x, 3 ) =qe*Light Speed* (-Q ( 1 , 1 ) ) *err(x,2); %vertical focussing force 

F (x, 4 ) =qe*Light Speed* (-Q ( 1 , 2 ) ) *err(x,2); % vertical defocussing force 

end 
end 



SHM Transfer Matrix (Ml.m): 

%SHM Matrix function 
function transfer = Ml (wb, T) 
transfer = [cos (wb*T) l/wb*sin (wb*T) 

-wb*sin (wb*T) cos (wb*T) ] ; 

end 
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2. Matlab program (ex4.m) for producing results in the form of: 




Main Program: 

Clear all 

% to simulate the effect of 84 kicks per revolution 

% on an electron in EMMA assuming uniform focusing strength 

% to produce graphs of maximum amplitude through N turns with relation to 

% increasing sigma values for the distribution of quadrupole alignment 

% errors (500 unique error patterns) . 

o, _ 
o 

NO = [100, 1000]; % no. of turns 

n=numel(N0); %number of elements in the array NO 

global Amax 

Amax = zeros (500,n) ; 

for turnset = l:n 

N = NO (turnset) ; 

ex4Nturns (turnset, N) ; 

end 

o, _ 
o 

%Graphing 

sigma = [1:500]' . * ( ( le-6) /5 ) ; %Sigma from 0.2 to 100 micrometers 
for b=l : n 
if b==l 

scatter ( sigma, Amax ( 1 : 500 , b) , 10, 'b') 
axis ( [0 le-4 0 0 . 025] ) 
xlabel ( ' 1 \sigma (m) ') 



66 



Chris Edmonds 



A Small Accelerator for Treating Cancer 



ylabel('SQRT (2J_{ z,max} ) [m A {l/2}] ') 
hold on 

p = polyfit (sigma, Amax ( 1 : 50 0 , b) , 1 ) ; %fits line of best fit to data 
lobf = p(l) . * sigma +p(2); 

plot (sigma, lobf) %plots line of best fit to data 
end 

if b==2 

scatter ( sigma, Amax ( 1 : 500 , b) , 20, ' rx ' ) 

p = polyfit (sigma, Amax ( 1 : 50 0 , b) , 1 ) ; %fits line of best fit to data 
lobf = p(l) . * sigma +p(2); 

plot (sigma, lobf, 'Color', 'red') %plots line of best fit to data 
end 
end 

hold off 



Program for Calculating the Maximum Amplitude Given a Rate of Acceleration and Standard 
Deviation of Alignment Errors: 



function ex4Nturns (turnset, N) 

global Q qe LightSpeed Amax 

CT = 42*N; %number of cells traversed 

L = 16.57; % (m) Circumference 

%Quadrupole details ( Column 1 is focussing and 2 is defocussing, 
%row 1 is field gradient and 2 is lattice length) 
lqf = 58.782e-3; %length of horizontally focusing quad 
lqd = 75.699e-3; %length of horizontally defocusing quad 
Q = [6.7, -4.9 
lqf, lqd] ; 

qe = 1.602e-19; ^electron charge 

LightSpeed = 2.998e8; % (m/s) speed of light 
^Assumes straight paths! 

TO = L/LightSpeed; %time taken for one revolution 

TDF = 0 . 274571/LightSpeed; %time taken to travel from defocussing to 
focussing quad 

TFD = 0 . 120023/LightSpeed; %time taken to travel from focussing to 
defocussing quad 

o, _ 
o 

% set matrices of before (1) and after (2) position and velocity 
xl = zeros (1, CT) ; % 

vl = zeros (1, CT) ; % xl (m) displacement & vl (m/s) velocity after one cell 
Al = zeros (1, CT) ; 
betal = zeros (1, CT) ; 
tunegraph = zeros (1,CT); 

o, _ 
o 

^Acceleration parameters 
%Energy range 

Estart = 10.75; %Injection energy MeV 
Estop = 20.75; ^Extraction energy MeV 

deltaE = (Estop - Estart) /CT; %energy increase per cell 
for nxsigmamin = 1:500 

o, _ 
o 

%dist ribut ion of misalignments and associated forces 
mu=0; %mean 

sigma= nxsigmamin* ( 1E-6 ) /5 ; %standard deviation 

seed = nxsigmamin A 2 ; %to provide repeatable values for testing 
F = gauerr (mu, sigma, seed); 
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o, _ 
o 

% Particle Start parameters 

tunemax = ( 0 . 054 9*Estart A 2 ) - (2 . 4572*Estart ) + 32.4479; 
freq = tunemax / TO; 
wb = 2*pi*freq; 

betal(l) = L/ ( 2 *pi*tunemax) ; % beta start value 

xl(l) = 0; % (m) start displacement 

vl(l) = 0; % (m/s) start velocity 

Al(l) = sqrt(xl(l) A 2 + vl(l) A 2 / wb A 2); 

%count ing 

c=l; 

for cl = 1 : (CT-1) , 

energy = Estart + cl*deltaE; 

tune = (0 . 0549*energy A 2) - (2 . 4572*energy ) + 32.4479; 
tunegraph (cl) = tune; 
freq = tune / TO; 
wb = 2*pi*freq; 

dvF = (F (c, 3) * LightSpeed * Q (2, 1) ) / (energy*le6*qe) ; ^transverse 
velocity kick at focussing quad 

dvD = (F(c,4) * LightSpeed * Q (2, 2) ) / (energy*le6*qe) ; ^transverse 
velocity kick at defocussing quad 

vectorl = [xl(cl) 

vl (cl) ] ; 

M = Ml (wb, TDF) ; 

vector2 = M * vectorl; 
following focussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + 
to vector 

M = Ml (wb, TFD) ; 

vector2 = M * vector2; 
defocussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + 
quad to vector 

xl(cl+l) = vector2(l); 
vl(cl+l) = vector2(2); 

betal(cl + l) = L/ (2^pi^tune) ; 

Al(cl + 1) = sqrt ( (xl (cl + 1) ) A 2 + ( vl ( ( cl + 1 ) ) A 2 / wb A 2)); 
c = c+1; 
if c>42 
c=l; 

end 
end 

% Maximum amplitude for nxsigmamin errors over N turns 
Al = Al ./ sqrt(betal); 
Amax (nxsigmamin, turnset) = max(Al); 
end 
end 



%travel from previous defocussing quad to 
dvF; %adds kick from horizontally focussing quad 

%travel from previous focussing quad to following 
dvD; %adds kick from horizontally defocussing 
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Random Cell Alignment Generator (gauerr.m): 

Refer to gauerr.m in the first programming section. 
SHM Transfer Matrix (Ml.m): 

Refer to Ml.m in first programming section. 
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3. Matlab program (ex5.m) producing a graph of the rate of change of maximum amplitude 
with respect to sigma Vs the rate of acceleration: 



700 




0 — 1 1 1 1 1 1 1 L 1 

0 2 4 6 8 10 12 14 16 18 

1/sqrt(dQ/dT) 

Main Program: 

clear all 

% to simulate the effect of 84 kicks per revolution 

% on an electron in EMMA 

% assuming uniform focusing strength 

% to produce graphs of the rate of change of maximum amplitude with respect 
% to sigma Vs the rate of acceleration. 

o, _ 
o 

NO = [100, 200,300, 400, 500,750, 1000,1500, 2000]; % no. of turns 

n=numel(N0); %number of elements in the array NO 

global Amax GradTune 

Amax = zeros (500, n); 

GradTune = zeros (l,n); 

for turnset = l:n 

N = NO (turnset) ; 

ex5Nturns (turnset, N) ; 

(turnset/n) *100 

end 

g, _ 
o 

sigma = [1:500]' . * ( (le-6) /5) ; 
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p=zeros (l,n) ; 
RMSD = zeros (1 , n) ; 

%fits line of best fit to data for each rate of acceleration and calculates 
%an error on the gradient of this line based upon the distribution of 
%points around it. 
for b=l : n 

fit4Nturns = polyfit (sigma, Amax ( 1 : 500 , b) , 1 ) ; 
p (l,b) =f it4Nturns (1) ; 

trend = f it 4Nturns ( 1 ) . *sigma + f it 4Nturns (2 ) ; 

Errl = 0; 

for errc=l : 500 

Errl = Errl + (Amax (errc , b) - trend(errc, 1) ) A 2; 
end 

RMSD(l,b) = sqrt (Errl/498) ; 
sumsigsq = sum ( sigma . A 2 ) ; 
sqsumsig = sum ( sigma) A 2 ; 

errgrad (1 , b) =RMSD (1 , b) * sqrt (500/ ( (500 * sumsigsq) -sqsumsig) ) ; 

end 

d = 1 . / sqrt (GradTune ) ; 

g, _ 
o 

%Least square fit 
Nu = numel (d) ; 
meand = sum(d)/N; 
meandsq= sum(d. A 2)/N; 

o, _ 
o 

%Gradient and constant for line in the form y=mx+c 

S = sum ( 1 ./ (errgrad . A 2 )) ; 

Sxy = sum ( (d . *p) ./ (errgrad . A 2 )) ; 

Sx = sum (d ./ (errgrad . A 2 )) ; 

Sy = sum (p ./ (errgrad . A 2 )) ; 

Sxx = sum ( (d . A 2 )./ (errgrad . A 2 )) ; 

delta = (S*Sxx) - (Sx A 2) ; 

m = ( (S*Sxy) - (Sx*Sy) ) /delta 

c = ( (Sxx*Sy) - (Sx*Sxy) ) /delta; 

o, _ 
o 

yd = zeros (l,Nu+l) ; 
yd (2 :Nu+l) =d (1 :Nu) ; 
y = m. ^yd + c; 

o, _ 
o 

%Error on gradient 
sigma = sqrt (S/delta) 

%Chi sq 

Chisq = sum (( (p-y ( 2 : Nu+1 ))./ (errgrad) ). A 2 ) 
ChisqNDF = Chisq/7 

o, _ 
o 

%Graphs 

errorbar (d, p, errgrad, 'bx') 
axis ( [0 18 0 900] ) 
xlabel ( ' 1/sqrt (dQ/dT) ' ) 

y label ( 1 \DeltaSqrt (2 J_{ z,max} ) / \Delta\ sigma [m A {-1/2 } ] ' ) 
hold on 
plot (yd, y ) 

o, _ 
o 
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Function for Calculating the Maximum Amplitude Given a Rate of Acceleration and Standard 



Deviation of Alignment Errors: 



function ex5Nturns (turnset, N) 
global Q qe LightSpeed Amax GradTune 
CT = 42*N; %number of cells traversed 
L = 16.57; % (m) Circumference 

%Quadrupole details ( Column 1 is focussing and 2 is defocussing, 
%row 1 is field gradient and 2 is lattice length) 
lqf = 58.782e-3; %length of horizontally focusing quad 
lqd = 75.699e-3; %length of horizontally defocusing quad 
Q = [6.7, -4.9 
lqf, lqd] ; 

qe = 1.602e-19; ^electron charge 

LightSpeed = 2.998e8; % (m/s) speed of light 
^Assumes straight paths! 

TO = L/LightSpeed; %time taken for one revolution 

TDF = 0 . 274571/LightSpeed; %time taken to travel from defocussing to 
focussing quad 

TFD = 0 . 120023/LightSpeed; %time taken to travel from focussing to 
defocussing quad 

o, _ 
o 

% set matrices of before (1) and after (2) position and velocity 
xl = zeros (1, CT) ; % 

vl = zeros (1, CT) ; % xl (m) displacement & vl (m/s) velocity after one cell 
Al = zeros (1, CT) ; 
betal = zeros (1, CT) ; 
tunegraph = zeros (1,CT); 

o, _ 
o 

^Acceleration parameters 
%Energy range 

Estart = 10.75; %Injection energy MeV 
Estop = 20.75; ^Extraction energy MeV 

deltaE = (Estop - Estart) /CT; %energy increase per cell 
for nxsigmamin = 1:500 

o, _ 
o 

%dist ribut ion of misalignments and associated forces 
mu=0; %mean 

sigma= nxsigmamin* ( 1E-6 ) /5 ; %standard deviation 

seed = nxsigmamin A 2 ; %to provide repeatable values for testing 
F = gauerr (mu, sigma, seed); 

o, _ 
o 

% Particle Start parameters 

tunemax = ( 0 . 054 9*Estart A 2 ) - (2 . 4572*Estart ) + 32.4479; 
tunemin = ( 0 . 054 9*Estop A 2 ) - (2 . 4572*Estop) + 32.4479; 
freq = tunemax / TO; 
wb = 2*pi*freq; 

GradTune ( 1 , turnset ) = (tunemax - tunemin) /N; %change in tune per turn, 
linear approximation. 

betal (1) = L/ ( 2 *pi*tunemax) ; % beta start value 

xl(l) = 0; % (m) start displacement 

vl(l) = 0; % (m/s) start velocity 

Al(l) = sqrt(xl(l) A 2 + vl(l) A 2 / wb A 2); 

%count ing 

c=l; 

for cl = 1 : (CT-1) , 

energy = Estart + cl*deltaE; 
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tune = (0 . 0549*energy A 2) - (2 . 4572*energy ) + 32.4479; 
tunegraph (cl) = tune; 
freq = tune / TO; 
wb = 2*pi*freq; 

dvF = (F (c, 3) * LightSpeed * Q (2, 1) ) / (energy*le6*qe) ; ^transverse 
velocity kick at focussing quad 

dvD = (F(c,4) * LightSpeed * Q (2, 2) ) / (energy*le6*qe) ; ^transverse 
velocity kick at defocussing quad 

vectorl = [xl(cl) 

vl (cl) ] ; 

M = Ml (wb, TDF) ; 

vector2 = M * vectorl; 
following focussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + 
to vector 

M = Ml (wb, TFD) ; 

vector2 = M * vector2; 
defocussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + 
quad to vector 

xl(cl+l) = vector2(l); 
vl(cl+l) = vector2(2); 

Al(cl+1) = sqrt (xl (cl+1) A 2 + (vl(cl+l) A 2 / wb A 2)); 
betal(cl+l) = L/ (2*pi*tune) ; 
c = c+1; 
if c>42 
c=l; 

end 
end 

% Maximum amplitude for nxsigmamin errors over N turns 
Al = Al . /sqrt (betal) ; 
Amax (nxsigmamin, turnset) = max(Al); 
end 
end 

Random Cell Alignment Generator (gauerr.m): 

Refer to gauerr.m in the first programming section. 
SHM Transfer Matrix (Ml.m): 

Refer to Ml.m in first programming section. 



%travel from previous defocussing quad to 
dvF; %adds kick from horizontally focussing quad 

%travel from previous focussing quad to following 
dvD; %adds kick from horizontally defocussing 
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4. Matlab program (ex6maxamp4startamp.m) for calculating the physical aperture: 
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Main Program 



clear all 



o. 
o 



o, 
o 



o. 
o 



o_ 
o 



o_ 
o 



to simulate the effect of 84 kicks per revolution 
on an electron in EMMA assuming uniform focusing strength 
with the aim of estimating the physical aperture for 100 turns 
and 10 0 0 turns . 



NO = [100, 1000]; % no. of turns 

n=numel (NO) ; %number of elements in the array NO 

global DA 

DA = zeros (500, n) ; 

DAplot = zeros (33, n) ; 

sigmaA= zeros (33) ; 

for turnset = l:n 

N = NO (turnset) ; 

ex6routine (turnset , N) ; 

end 



sigma = [1:500]' . * ( ( le-6 ) /5 ) ; 
for b=l : n 



o_ 
o 



Amid = -7; 
for A = 1 : 33 

Amid = Amid+15; 
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sigmaA(A) = sigma (Amid) ; 
avg = 0 ; 

for range = -7:7 
f =Amid+ range ; 
avg = avg + DA(f,b); 

end 

avg = avg/ 15; 
DAplot (A, b) = avg; 

end 

o, _ 
o 

if b==l 

% scatter (sigma, DA (1:500, b), 20, ' bx ' ) 
sigl = sigmaA ( 1 : 33 ) ; 
Dal = DAplot (1:33,1); 

semilogy ( sigl , Dal , 'blue', ' LineWidth ' , 3 ) 
xlabel ( ' 1 \sigma (m) ') 

ylabel (' Dynamic Aperture (\pi m rad) ') 
axis([0 le-4 10e-7 le-4]) 
hold on 
end 

if b==2 

% scatter (sigma, DA (1:500, b), 20, ' rx ' ) 
sig = sigmaA ( 1 : 1 7 ) ; 
Da = DAplot (1:17,2); 

semilogy ( sig, Da, 'red', ' LineWidth ', 3 ) 
end 
end 

hold off 

Function for Calculating the Maximum Start Amplitude of a Particle which is Accelerated to 
the Extraction Energy Without Incidence Upon Accelerator Wall Given A Particular Standard 
Deviation of Error Distribution and Rate of Acceleration 

function ex6rout ine ( turnset , N) 
global Q qe LightSpeed DA 

o, _ 
o 

^Accelerator Constants 

CT = 42*N; %number of cells traversed 
L = 16.57; % (m) Circumference 

%Quadrupole details ( Column 1 is focussing and 2 is defocussing, 

%row 1 is field gradient and 2 is lattice length) 

lqf = 58.782e-3; %length of horizontally focusing quad 

lqd = 75.699e-3; %length of horizontally defocusing quad 

Q = [6.7, -4.9 

lqf, lqd] ; %Quadrupole gradients and lengths 

qe = 1.602e-19; ^electron charge 

LightSpeed = 2.998e8; % (m/s) speed of light 
^Assumes straight paths! 

TO = L/LightSpeed; %time taken for one revolution 

TDF = 0 . 274571/LightSpeed; %time taken to travel from defocussing to 
focussing quad 

TFD = 0 . 120023/LightSpeed; %time taken to travel from focussing to 
defocussing quad 
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o. 
o 



% set matrices of position and velocity 
xl = zeros (1, CT) ; % 

vl = zeros (1, CT) ; % xl (m) displacement & vl (m/s) 



velocity after one cell 



o_ 
o 



^Acceleration parameters 
%Energy range 

Estart = 10.75; %Injection energy MeV 
Estop = 20.75; ^Extraction energy MeV 

deltaE = (Estop - Estart) /CT; %energy increase per cell 
SA=[0: 0.1:11] . *le-3; %Start amplitudes 

tunestart = (0.0549* (Estart A 2) ) - (2 . 4 5 72*Estart) +32.447 9; 
beta = L/ ( 2 *pi*tunestart ) ; 
root2 J=8 . 9e-3/sqrt (beta) ; 



o_ 
o 



%dist ribut ion of misalignments and associated forces 

mu=0; %mean 

for nxsigmamin = 1:500 

nxsigmamin 

sigma= nxsigmamin* ( 1E-6 ) /5 ; %standard deviation 

seed=nxsigmamin A 2 ; 

F = gauerr (mu, sigma, seed); 

k=0; 

start amp = 0; 
while k == 0 

startamp = startamp +1; ^increases start amplitude from zero to half 
aperture 

o, _ 
o 

% Particle Start parameters 

xl(l) = SA ( startamp) ; % (m) start displacement 
vl(l) = 0; % (m/s) start velocity 

%count ing 

c=l; %individual cell number (1-42) 

o, _ 
o 



%Calculate maximum amplitude for trial start conditions 
Al = 0; 
JA1=0; 
cl = 0; 

while JAKroot2J && cl<=(CT-l) %while amplitude is less than 
aperture size and acceleration is not complete 
cl=cl+l; %counts total number of cells 
energy = Estart + cl*deltaE; 

tune = (0.0549* (energy A 2) ) - (2 . 4 5 72*energy) +32.447 9; 
freq = tune / TO; 
wb = 2*pi*freq; 

dvF = (F(c,3) * LightSpeed * Q (2 , 1 ) ) / (energy*le6*qe) ; ^transverse 
velocity kick at focussing quad 

dvD = (F (c, 4) * LightSpeed * Q (2 , 2 ) ) / (energy*le6*qe) ; ^transverse 
velocity kick at defocussing quad 



vect or 1 



[xl (cl) 
vl (cl) ] ; 



M = Ml (wb, TDF) ; 

vector2 = M * vectorl; %travel from previous defocussing quad to 
following focussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + dvF; %adds kick from horizontally focussing 
quad to vector 
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M = Ml (wb, TFD) ; 

vector2 = M * vector2; %travel from previous focussing quad to 
following defocussing quad 

vector2 ( 2 ) =vector2 ( 2 ) + dvD; %adds kick from horizontally 
defocussing quad to vector 



xl(cl+l) = vector2(l); 
vl(cl+l) = vector2(2); 

Al = sqrt (xl (cl+1) A 2 + (vl(cl+l) A 2 / wb A 2)); 
betal = L/ (2*pi*tune) ; 
JA1 = Al/sqrt (betal) ; 
c = c+1; 

if c>42 

c=l; 
end 

end 

g, _ 
o 



DA (nxsigmamin, turnset) = ( ( SA ( startamp) A 2 ) ) / (pi* (beta) ) ; 
if JA1 > root2J 
k=l; 

end 



end 
end 
end 



Random Cell Alignment Generator (gauerr.m): 



Refer to gauerr.m in the first programming section. 



SHM Transfer Matrix (Ml.m): 



Refer to Ml.m in first programming section. 
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